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I.  Introduction 


The  process  of  mixing  of  gaseous  fuel  and  oxidant  is  one  of  fundamental  subjects  of 
combustion  and  propulsion  research.  It  is  basically  clear  that  the  better  (faster)  mixing  the  higher 
combustion  efficiency  can  be  reached  under  similar  conditions.  The  mixing  can  be  critically 
important  for  short  residence  time  combustion  systems  such  as  scramjet. 

The  theory  of  combustion  in  so-called  premixed  systems  is  developed  for  a  long  time  [1-3]. 
The  key  question  in  such  a  system  is  the  flame  propagation. 

Two  different  modes  of  the  flame  propagation  are  usually  considered:  laminar  and  turbulent. 
The  laminar  process  of  combustion  is  defined  by  balance  of  heat  release  in  the  flame  front  due  to 
chemical  reactions  between  fuel  and  oxidant  and  energy  transfer  by  heat  conductivity. 

In  the  case  of  turbulent  combustion  the  effective  flame  propagation  velocity  increases 
significantly  due  to  important  role  of  convection  in  energy  transport  from  the  heat  release  region  to 
the  undisturbed  flammable  mixture.  In  fact  the  turbulence  in  the  combustion  zone  can  be  created 
due  to  intensive  heat  release  (heating  resulting  in  non-steady  state  expansion  flow  of  gaseous 
products  and  fresh  mixture)  or/and  it  can  be  convected  by  the  flow  to  the  combustion  zone  as  well 
[3]. 

In  contemporary  combustion  systems  such  as  combustion  chamber  of  (turbo)jets,  ramjet, 
scramjet  combustion  starts  in  only  partially  mixed  flow  where  the  fuel/oxidant  ratio  varies 
significantly  in  working  volume.  In  this  particular  but  very  important  practical  case  the  combustion 
and  mixing  can  develop  simultaneously  affecting  each  other  [4-6]  (see  also  [7]  with  many  example 
of  particular  recently  performed  studies  in  the  field).  The  combustion  as  an  elementary  chemical 
reaction  between  fuel  and  oxidant  molecules  needs  the  direct  contact  (collision)  of  these  molecules. 
Thus,  the  combustion  can  occur  if  in  a  particular  location  the  collision  between  fuel  and  oxidant 
molecule  can  take  place.  This  «molecular»  level  of  mixing  is  provided  by  conventional  molecular 
diffusion.  However,  the  molecular  diffusion  is  generally  speaking  very  slow  process  in  scale  of 
available  mixing  time,  which  is  limited,  by  residence  time  in  combustion  zone.  To  enhance  the 
mixing  process  the  other  «large  scale»  mechanism  of  mixing  is  introduced.  The  main  idea  is  to 
increase  the  contact  surface  between  fuel  and  oxidant.  The  simplest  well-known  method  is  rotating 
mixing  widely  used  in  mankind  practice  for  whole  human  history.  The  turbulence  motion  is 
basically  a  vortex  motion  and  one  can  consider  the  turbulence  enhancement  of  mixing  as  a 
combined  effect  of  different  scale  rotating  mixing  cells.  This,  kinematics ,  stage  of  mixing  increases 
the  all  inter-component  contact  surfaces.  (The  theory  of  kinematics  mixing  based  on  a  Lagrangian’ 
non-turbulence  approach  is  recently  presented  in  a  textbook  by  J.M.Ottino  [8].  The  kinematics 
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treatment  of  the  turbulence  for  the  mixing  study  is  independently  proposed  and  developed  in  our 
recent  papers  [9- 12].) 

Simultaneously  the  conventional  molecular  diffusion  across  the  contact  surface  develops 
providing  the  mutual  penetration  of  fuel  into  the  oxidant  fluid  domain  and,  the  oxidant,  to  the  fuel 
fluid  domain.  As  a  result  at  the  contact  fuel  and  oxidant  domain  a  intermediate  diffusion  region  is 
formed  where  the  fuel-to-oxidant  ratio  varies  between  its  limits  0  and  «>.  The  characteristics 
thickness  of  such  a  layer  is  estimated  by 
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where  D  and  t  are  diffusion  coefficient  and  time  of  the  contact,  correspondingly. 

Thus,  increasing  contact  surface  and  diffusions  layer  forms  a  reacting  volume  where  the 
combustion  can  occur. 

The  concept  of  the  reacting  volume  has  been  formulated  by  us  in  series  of  papers  dealing 
with  the  evolution  of  plasma  clots  in  a  channel  flow  [9-11],  and  more  recently  to  treat  the  current 
break-off  in  a  gaseous  current  breaker  device  [12].  In  these  studies  the  reacting  volume  evolution 
was  described  with  taking  into  account  of  the  turbulence  flow  field.  Furthermore,  the  turbulence  is 
considered  as  flow  characteristics  modified  by  electromotive  body  force  jxB.  The  important  feature 
of  the  body  force  is  its  vortex  nature. 

This  concept  was  rather  fruitful  and,  for  this  reason,  it  was  proposed  to  check  a  similar 
approach  to  analysis  of  non-premixed  combustion  in  supersonic  flows  with  enhanced  turbulence  by 
electromotive  body  force  jxB. 

The  objective  of  the  present  work  is  (1)  to  check  the  reacting  volume  concept  for 
combustion  study  applications;  (2)  to  adopt  the  numerical  models  developed  for  gas-plasma  flows 
analysis  for  the  combustion  conditions;  (3)  to  develop  the  numerical  model  of  combustion  under 
high-speed  flow  conditions  with  MHD  interaction  assisted  fuel-oxidant  mixing. 
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II.  Kinematics  of  mixing  in  fluids 


Presently,  the  study  of  fluid  mixing  has  very  little  scientific  basis;  processes  and  phenomena 
are  analyzed  on  a  case-by-case  basis  without  any  attempt  to  discover  generality.  For  example,  the 
analysis  of  mixing  and  ‘stirring’  of  contaminants  and  tracers  in  two-dimensional  geophysical  flows 
such  as  in  oceans;  the  mixing  in  shear  flows  and  wakes  relevant  to  aeronautics  and  combustion;  the 
mixing  of  fluids  under  Stokes’s  regime  generally  encountered  in  ‘blending’  of  viscous  liquids  such 
as  polymers;  and  the  mixing  of  diffusing  and  reacting  fluids  encountered  in  various  types  of 
chemical  reactors  share  little  in  common  with  each  other,  except  possibly  the  nearly  universal 
recognition  among  researchers  that  they  are  very  difficult  problems. 

There  are,  however,  real  similarities  among  the  various  problems  and  the  possible  benefit 
from  overall  attack  on  the  problem  of  mixing  using  a  general  viewpoint  are  substantial.  The 
contribution  of  the  present  very  particular  efforts  can  be  useful  in  the  field  of  constructing  of 
practical  evaluation  procedures,  which  can  provide  more  direct  references  to  the  experimental  data, 
and,  in  turn,  development  advanced  mixing  technique. 

The  point  of  view  adopted  in  this  study  based  on  two  main  ideas. 

The  first  one  seems  to  be  fully  similar  to  the  formulated  by  J.Ottino  in  his  book  [8]  is  that 
from  kinematics  viewpoint  fluid  mixing  is  the  efficient  stretching  and  folding  of  material  lines  and 
surfaces.  Such  a  problem  corresponds  to  the  solution  of  the  dynamical  systems 

dx/dt  =  v(x,t) 

with  x  =  X  at  time  f=0. 

The  solution  of  the  above  system  is  called  the  flow  motion  in  the  approach  of  the 

book[8].  It  is  explained  there  that  even  in  the  simplest  cases  this  “solution”  might  be  extremely  hard 
to  obtain.  The  kinematical  foundations  lie  in  an  understanding  of  the  point  transformation  • 

In  the  frame  of  such  an  approach  the  following  three  sub-programs  are  considered: 

(1)  Mixing  of  s  single  fluid  or  similar  fluids. 

The  basic  objective  here  is  to  compute  the  length  (of  area)  corresponding  to  a  set  of  initial 
conditions.  As  it  can  be  shown  (see,  for  example,  [8])  only  in  few  cases  it  can  be  done  exactly  and 
in  most  of  these  the  length  is  mild.  The  best  achievable  mixing  corresponds  to  exponential 
stretching  nearly  everywhere  and  occurs  in  some  regions  of  chaotic  flows. 
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It  should  be  noted  also  that  the  latter  case  was  considered  in  our  papers  [9-13]  as  main  mechanism 
of  the  mixing  in  turbulent  flows  with  the  reference  to  Batchelor’s  hypothesis  of  exponential 
stretching  law. 

However,  under  these  conditions  the  (exact)  calculation  of  the  length  and  location  of  lines 
and  areas  is  hopelessly  complicated.  As  it  is  shown  in  [8]  even  extremely  simplified  flows  might  be 
inherently  chaotic  and  a  complete  characterization  is  not  possible.  For  example,  from  a  dynamical 
systems  viewpoint  if  the  system  possesses  horseshoes  infinitely  many  periodic  points  occurred  all 
of  which  cannot  be  calculated  precisely.  Fortunately,  as  far  as  mixing  is  concerned  the  low  periodic 
events  are  of  interest,  since  the  quick  mixing  is  desirable.  Nevertheless,  there  is  always  the  intrinsic 
limitation  of  being  unable  to  calculate  precise  information  (most  practical  problems  involve 
stretching  of  order  104  or  higher)  such  as  length  stretch  and  location  of  material  surface.  This  is  true 
even  though  none  of  the  flows  discussed  in  [8]  is  turbulent  in  an  Eulerian  sense.  Rather,  the 
previous  finding  should  be  used  to  establish  the  limits  of  what  might  constitute  reasonable  answer 
in  more  complicated  flows  (real  turbulent  flows  come  immediately  to  mind).  The  problem  here  how 
to  best  characterize  the  mixing,  knowing  beforehand  that  a  complete  characterization  is  impossible 
[8].  Furthermore,  the  practical  goal  is  not  only  to  characterize  the  mixing  but  also  rather  to  control 
the  mixing!  So  the  practical  approach  developed  in  this  project  involves  the  second  idea  of  reacting 
volume  which  will  be  described  later. 

(2)  Within  a  framework  of  a  family  flows  xs=<5?sX(X),  s=  each  of  motions  is  assumed  to 

be  topological  (see  Section  2.3  of  [8]);  mixing  of  similar  diffusing  and  reacting  fluids. 

This  case  corresponds  to  the  case  of  mixing  of  two  streams,  composed  of  possible  several 
species  that  are  Theologically  identical;  i.e.  they  have  the  same  density,  viscosity,  etc.,  and  have  no 
interfacial  tension.  Concurrently  with  the  mechanical  mixing  there  is  mass  diffusion,  and,  possible, 
chemical  reaction.  However,  for  simplicity,  it  is  typically  assumed  that  neither  the  diffusion  nor  the 
reaction  affects  the  fluid  motion.  This  case  is  discussed  in  detailed  in  [8]  as  a  lamellar  structures 
case. 

Again,  the  implementation  of  such  a  kinematics  approach  to  the  particular  problems  considered  in 
the  project  needs  an  accounting  in  some  way  of  back  affects  of  the  (( non-linear  ”  mixing 
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consequences  on  the  flow  field  (fluid  motion).  Furthermore,  some  kind  of  controllable  back  affects 
can  be  utilized  in  principle  as  an  intensification  of  mixing  -  main  goal  of  he  present  efforts. 

(3)  Mixing  of  different  fluids;  case  in  which  the  motions  are  non-topological,  i.e.,  there  is 
breakup  and/or  fusion  of  material  elements. 

In  this  case,  the  mixing  of  two  or  more  fluids  leads  to  breakup  and  coalescence  of  material 
regions.  This  problem  is  complicated  and  only  few  special  cases  belonging  to  this  category  are 
discussed  in  [8]. 

In  Fig.l  through  5  several  illustrations  of  the  kinematics  mixing  are  presented  to  be 
compared  later  with  our  simulation  results  obtained  in  course  of  the  gas-plasma  flow  studies. 


Fig.l.  Basic  processes  occurring  during  mixing:  (a)  corresponds  to  the  case  of  two  similar  fluids 
with  negligible  interfacial  tension  and  negligible  inter-diffusion:  an  initially  designated 
material  region  stretches  and  folds  by  the  action  of  a  flow;  (b)  corresponds  to  a  blob  diffusing 
in  the  fluid:  in  this  case  the  boundaries  become  diffuse  and  the  extent  of  the  mixing  is  given 
by  level  curves  of  concentrations  (a  profile  normal  to  the  striations  is  shown  at  the  right);  in 
(c)  the  blob  breaks  due  to  interfacial  tension  forces,  producing  smaller  fragments  which  might 
in  turn  stretch  and  break  producing  smaller  fragments.  Case  (b)  is  an  excellent  approximation 
to  (a)  if  diffusion  is  small  during  the  time  of  the  stretching  and  folding.  In  (a)  the  blob  is 
passive,  in  (c)  the  blob  is  active. 
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Fig.2.  Growth  of  a  “material  line”  composed  of  small  hydrogen  bubbles  produced  by  a  platinum 
wire  stretched  across  an  decaving  turbulent  flow  behind  a  grid  placed  at  the  extreme  left.  The 
Reynolds  number  based  on  the  wire  diameter  is  1.360.  (Reproduced  with  permission  from 
Corrsin  and  Karweit  (1969).) 


Fig.3.  Mixing  in  a  creeping  flow.  The  figure  shows  the  deformation  of  a  material  region  in  a  journal 
bearing  flow  when  it  operates  in  a  time-periodic  fashion  (experiment  from  Chaiken  et.al. 
(1986)). 


Fig.4.  Concentration  of  a  turbulent  round  jet  into  water  at  Reynolds  number  2.300,  measured  by 
laser  induced  fluorescence;  the  cut  is  along  a  plane  including  the  axis  of  a  symmetry  of  the  jet. 
(Reproduced  with  permission  from  Dimotakis,  Miake-Lye,  and  Papantonionu  (1983).) 
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Fig.5.  Deformation  of  a  tracer  in  a  numerical  experiment  of  motion  in  the  Earth  mantle.  The  sides  of 
the  rectangle  are  insulating  but  the  bottom  is  subjected  to  a  constant  heat  flux  while  the 
temperature  of  the  top  surface  is  kept  constant.  The  motion  is  produced  by  buoyancy  and 
interval  heating  elects  (the  fluid  is  heated  half  from  below  and  half  from  within).  The 
Rayleigh  number  is  1.4xl06,  the  time  scale  of  the  numerical  simulation  corresponds  to  155 
Myear,  and  the  thickness  off  the  layer  is  700km.  An  instantaneous  picture  of  the  streamlines 
reveals  five  cells.  (Reproduced  with  permission  from  Hoffman  and  McKenzie  (1985).) 


III.  Gas-Plasma  Flow  Examples 


On  a  Current  Carrying  Plasma  Clot  Evolution  in  Turbulent  Flow  [11-12] 


Current  carrying  clot  created  in  high  velocity  flow  in  an  MHD  channel  is  influenced  by 
many  factors.  Energy  balance  of  the  clot  is  mainly  defined  by  Joule  dissipation  inside  the  clot  and 
energy  losses  trough  its  borders.  The  analyses  of  the  current  carrying  clot  energy  balance  in 
turbulent  flow  performed  in  our  papers  [11-12]  has  shown  that  the  short  lifetime  features  of  the 
clots  are  strongly  important  for  the  clot  evolution.  It  was  found  also  that  the  clot  temperature 
evolution  depends  on  the  initial  shape  and  turbulence  parameters  of  the  flow.  Three  different  type  of 
clot  temperature  evolution  were  described.  The  electrode  wall  losses  and  the  non-potential  body 
force  jxB  are  found  to  be  very  important  as  well. 

1.  Turbulent  deformation  of  plasma  clot. 

The  evolution  of  the  initially  compact  volume  of  fluid  Vm  limited  by  surface  S  in  the  flow  of 
incompressible  fluid  with  uniform  isotropic  turbulence  has  been  discussed  in  Ref[9]. 

Namely,  it  was  shown  in  Ref.[9]  that  the  surface  growing  is  defined  by  the  relation 

S(0  =  Soexpi-*7*E(*)rf4  (1) 

lM  0  J 

In  this  way  surface  seems  to  depend  on  time  by  exponent,  and  this  is  in  accordance  with  Ref 
[14]  For  evaluation  of  factor  in  the  exponent  power  one  should  take  E(k)  from  Kolmogorov’s 
spectrum  (see,  for  example  Ref[15,16]: 

E(k)  =  t2/3k~5/3 

where  £  is  turbulent  energy  dissipation  £~v3//,  v  is  basic  scale  of  fluctuation  velocity,  /  is  basic  scale 
of  turbulence  [17].  Integrating  in  (1)  we  get: 
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S(t )  =  Sq  exp 


-■tz2^ 


J* 


-2/3 


dk 


(2) 


Taking  instead  of  kmm  following  to  Ref.[17]  (Re/Recrit )3/4//  (Re  is  Reynolds  number  of  the 
basic  motion  and  Recrjt  is  Reynolds  number  of  transmit  from  laminar  to  turbulent  motion),  instead  of 
km\n  1//  and  considered  that  (Re/Recrit)3/4»1  we  get  finally: 


e_o  v  (Re/Recrt7)  t  ,,, 

S  —  Sq  exp<  ~  r  wf 

It  is  obvious,  that  in  this  kind  formulas  relations  between  values  u,  v  and  /  are  define  with 
some  uncertainty.  For  this  reason  some  factor  of  order  of  1  would  be  included  in  the  power  of 
exponent  and  from  the  other  hand  factors  ~  1  in  pre-exponent  expression  can  be  canceled  in  this 
estimation. 

Thus,  the  relation  (2)  should  be  used  in  practice  in  the  form: 

S-SoexpK^^1'4' 

ul 


where  factor  A  should  be  chosen  by  comparing  with  experimental  data. 

It  is  easy  to  see,  that  the  turbulence  «accelerates»  any  processes,  concerned  with  fluxes  of 
mass,  momentum  and  energy  (heat  transfer,  diffusion  and  so  on)  between  volume  Vm  and  ambient 
fluid  turbulence  by  means  this  kind  of  the  surface  increasing. 

Note  also,  that  one  can  estimate  contribution  of  any  scale  of  turbulent  motion  in  increasing 
of  the  surface  to  take  in  formula  (3)  integration  limits  k\  and  ki,  with  Ak=k\-k2  of  order  &i: 


S(t)  -  Sq  exp- 


Furthermore  the  turbulent  diffusion  D\  can  be  also  evaluated  by  using  this  fact.  It  depends 
mainly  on  maximal  scale,  that  for  given  time  t  take  part  in  the  surface  deformation  (for  example 
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which  increase  initial  surface  in  e  times),  i.e.  D\  ~  A,]Vi,  where  Xi~l/Ari  and  v\  fluctuation  velocity 
for  scale  X-i.  For  sufficiently  large  time  t  we  get  well-known  formula:  D  ~  Iv,  because  of  £mj„~l//. 


2.  Energy  Effects  of  Turbulence 


In  this  section  evaluation  of  the  effect  of  turbulence  on  the  energy  state  inside  a  current 
carrying  clot  will  be  obtained  when  current  passed  through  it. 

A  compact  clot  in  incompressible  fluid  volume  Vm=Sxl]Xl2  with  length  11  approximately 
equal  to  cross  section  perimeter  13  is  considered.  It  is  assume  that  the  clot  created  in  the  flow  is 
uniform  at  the  initial  time  moment  t= 0  with  constant  temperature  Tq»Tc  where  Tc  is  ambient 
temperature.  The  electrical  power  is  supplied  to  the  clot  by  constant  electromotive  force  E  power 
source  schematically  shown  in  Fig.6. 


Fig.6.  Schematic  discharge  structure  in  channel  flow. 


During  the  process  of  motion  turbulence  results  in  the  clot  deformation  to  some  layers  with 
the  thickness  5 (/)  and  the  side  surface  of  clot  S(t).  Time  dependence  for  S{t)  from  (3): 

S(t)  =  S0  exp  {Clt} , 

where  £2  is  characteristic  frequency  of  surface  deformation  development  and  in  terms  of  Section  1 : 


1/4 


(4) 


13 


By  virtue  of  the  clot  volume  resistance  estimation  from  [11]  the  Joule  dissipation  value 


reads: 


q  =  oo(7^So=ao(T)E2exp{_Qt}  (5) 

o(T) 

In  respect  to  diffusion  the  whole  region  considered  can  be  divided  in  three  sub  regions.  The 
first  is  the  region  of  cold  surrounded  gas,  the  second  is  the  region  occupied  by  the  clot  fraction 
being  not  affected  by  diffusion,  and  the  third  is  the  mixing  region  with  width  Ip  ~  <J~Dt .  Width  of 
this  region  changes  in  time  according  to  said  above. 

In  this  way  the  equation  of  energy  is  the  following: 

P  Cp¥-  =  q-kTV2r  (6) 

Using  (5),  the  dependence  of  conductivity  on  temperature  from  [18]: 


a(T)  =  Co  exp-!  ml 


To' 


J_ 

T 


and  evaluating  VT  as  THd  where  T  is  the  clot  temperature,  the  equation  (6)  can  be  presented  in  the 
form 


ft 

— =  5exp  {m{  1  -1/0)  -x}-XT -j - fD , 

IfoCpQ 


(7) 


where 


_Jl,  8</d 

/d  jo,  5>lD ' 

Here  new  variables  Q  =  T/Tq  and  t  =  tQ,  and  new  designation 
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(8) 


B  __  QqEq 
P  CpQT0 


Following  [1 1,12]  one  obtains  for  the  small  values  of  t 

0(t)  =  — — — exp{(w5-l)x} +  1  (9) 

mB- 1 

The  dependence  of  the  temperature  0  upon  time  ris  shown  in  Fig.7. 


Fig.7.  Bifurcation  of  linerized  solution  of  energy  equation  at  the  initial  point. 


One  can  see  that  in  the  initial  moment  three  cases  are  available  that  are  controlled  by  the 
parameter  mB-l.  If  mB>  1  then  the  clot  grows  exponentially;  if  mB- 1  then  this  growth  is  linear,  and 
if  mB<  1  then  the  temperature  quickly  comes  to  the  saturation.  During  the  time  all  curves  come  to 
saturation  (Fig.8)  but  its  level  is  quite  different. 

One  can  note  also  that  if  clot  form  is  a  stretched  cylinder  the  geometrical  factor  is  to  be 
taken  into  account 

Considering  such  a  cylinder  as  a  number  of  elements  regarded  above  it  needs  to  take  instead 

of  (17): 


b'=b-i2/l 2 
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where  /  is  the  width  of  cylinder  and  L  is  its  length. 


Fig.  8.  Dimensionless  clot  temperature  0  vs  time  T. 


3.  Electrode  Effects 

In  order  to  estimate  an  electrode  wall  effect  on  the  energy  balance  due  to  resistive  layer  near 
«cold»  wall,  the  resistance  of  the  boundary  layer  /?/,/  is  introduced  in  the  expression  for  Joule  heat 
release  in  the  form: 


q  = 


Ro+Rbi 


Ro- 


(10) 


Than  the  energy  equation  (7)  is  transformed  to  the  form 
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30 

3x 


B  exp  {m(l  - 1/0)  - 1} 


1  +  Rbia°VW-  exp  {w(l  - 1/0)  -  T} 
s0 


2  ~^T~2 - ~fD 

^2  ihcpo. 


The  dependence  of  the  temperature  q  upon  time  t  calculated  with  this  equation  is  shown  in 

Fig.9. 


Fig.9. 


A  linearized  expression  similar  to  (9)  is  now  as  following 

m=B - explore}  -1}  t|  (U) 

co  (1  +  £)(1  +  %  exp  {cox} ) 


where 

0)  =  _  l  fi  _ 

(£+l)2  ’  so 

Saturation  with  level  0um,  depended  upon  can  be  seen  in  Fig.  1 0. 
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B  1 


Slim© 


1  + 


K)|  1+5 


1 ,  co(^)  <  0  ’ 


Fig.  10. 


i.e.  quickly  decrease  with  increasing  E,.  Thus  taking  into  account  the  wall  effect  results  in  limitation 
of  exponential  growth  of  the  temperature. 
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4.  MHD  Generation  of  Turbulence. 


The  model  described  above  regards  the  effect  of  turbulence  coming  into  flow  from  outside, 
but  if  there  is  some  external  force  turbulence  is  bom  in  process  of  movement  directly  by  non¬ 
potential  part  of  one.  Further  evaluations  permit  to  define  characteristic  scale  of  turbulence, 
fluctuation  velocity  and  accordingly  coefficient  of  turbulent  diffusion.  We  assume  for  simplify  that 
fluid  is  incompressible. 

Equation  for  vorticity 

^>=ro(vxrotv]+^I  (12) 

dt  p 

where  f  is  external  volume  force,  and  v  is  fluctuation  of  velocity. 

Multiplying  by  rot  v  on  the  left  (in  this  case  first  term  on  the  right  is  equal  zero)  integrate  in 
the  volume  for  getting  of  average  value 

((rotv)2)  =  -i{(rotv)2JF 
flfa(rotv)2Lr_  rfotvrotf 

V2{  *  J  V  P 

_  rfotv  rotrotf  dy  |  |-[vxrotf]„  ^ 

V  p  s  p 

Hear  the  formula  of  vector  analysis  div[axb]  =  (brota)-(arotb)  and  Gauss  theorem  has 
been  used. 

Assume,  that  f  is  stepwise  function,  that  is: 

ff,  insideoftheclot 
[0,  elsewhere 

we  get  that  first  term  on  the  right  is  equal  zero. 

In  this  way  by  order  value: 
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(13) 


d  <  (rot  v)2  >  _  Sv  8f 

dt  pvT 

where  S,  r  are  characteristic  bound  surface  and  density,  l  is  characteristic  length  of  force  jump,  8f  (= 
f  in  our  case)  is  jump  itself  and  v  is  characteristic  fluctuation  velocity. 

Equation  for  velocity  we  can  present  in  such  form: 


d(vo  +  v) 

dt 


=  -V 


f 

(  i'\ 

P0.I0 

p  V 

Ip 

2J 

[p+TJ 

-V(v0v)  +  - 
P 


(14) 


where  Vo,  po  are  velocity  and  pressure  of  fluid  movement  v  is  fluctuation  of  velocity  and  p  is 
fluctuation  of  pressure. 

Multiplying  by  v  on  the  left,  canceling  cross  products  between  basic  quantities  and 
fluctuation  quantities  we  integrate  for  getting  <V2>: 


1  1  3<  v2  > 

2  V  dt 


+ 


J 


-  JvV(vqv)^F  +  j—dV 
V  V  p 


First  two  terms  on  the  right  are  common  gas  dynamics  terms,  but  we  have  interest  only  in 
the  last  one,  which  is  concerned  with  volume  force  and  hence  we  remain  only  this  term. 

Thus: 


3<t2>.£  (15) 

dt  p 

Now  evaluate  turbulent  movement  scale. 

For  this  purpose,  divide  square  of  vorticity  by  square  of  velocity,  regarded  that  characteristic 
times  of  value  change  are  the  same. 
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(16) 


,2_JvSv£_lV 

p  IV p  s 

Average  fluctuation  velocity  we  get  from  (15),  left  side  of  which  is  turbulent  energy 
dissipation. 

In  this  manner:  _/v/r=e=v3//,  therefore 


(17) 


For  diffusion  coefficient: 


l3V3  f 2 
Z)  =  X.v  =  j  ~3_2- 

I  P2 


1 1/4 


(18) 


Returning  to  integrate  volume  and  integrate  surface  can  be  gotten: 
characteristic  size  of  the  clot)  S~L2,  thus: 


l}l<V<I?  (L  is 


12<12<IL 


V'2 


<v< 


VK7 


r  2V/4 

uLi 

P  J 


'Hi"'2 

V  P  ^ 


<D< 


fL5l5f2^'4 


v 


(19) 


Note  also,  that  in  above  formulas  dependents  on  the  time  is 


implicit,  resulted  from  /  ~  -jDt . 


5.  Energy  Effects  of  Turbulence  in  Presence  of  Magnetic  Field. 

Let’s  consider  again  the  problem  solved  in  Section  2  in  presence  of  magnetic  field  B 
perpendicular  to  electrodes. 
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Again  shall  search  the  decision  for  the  region  occupied  by  the  clot  fraction  being  not 
affected  by  diffusion. 

Then  energy  equation  is: 

where  Eo  —  external  electric  field  and  j  —  volume  average  current  density. 

From  section  2  follows,  that  it  is  possible  to  define  j  as: 

j  =  cQ(T)(E0-uB)eQt  (21) 

where  u  —  average  flow  velocity. 

In  this  work  Q  represents  surface  growth  due  to  external  turbulence  as  well  as  induced 
turbulence  which  is  bom  by  nonpotential  force  jxB. 

In  above  section  average  fluctuation  velocity  scale  v  and  turbulent  movement  scale  1  were 
estimated.  In  addition  integral  character  of  this  and  by  virtue  of  integrated  character  of  these 
parameters  we  can  canceled  from  ideal-step  distribution  of  force  amplitude.  Having  entered  an 
assumption,  that  the  dependence  Qt  —  characteristic  of  induced  turbulence  on  v  and  /  like 
dependence  for  external  can  be  received: 


where /  —  non-potential  force. 

By  virtue  of  independence  of  these  two  processes  of  surface  change  it  is  possible  to  present 
£2  as  a  sum: 

£2  —  £20  +  £2  y , 


where  £2o  —  characteristic  of  external  turbulence. 

For  simplicity  now  let's  assume,  that  u= 0  and  unites  (1)  and  (2)  in  a  system,  previously 
having  reduced  them  to  dimensionless  variable: 
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/  =  em(l-l/e)e(l+a/)T 

30  n. 


Here: 

J  = 

7  . 

Ij  ;  0  =  - 

27o  : 

a  = 

OqEqB 

_  GqEqB 

P«Qq 

pv,2A, 

P  = 

°o£'o  . 

Jl  _  =1;  0l  , 

l/=0  5 

Expressing  q  from  (22*)  through  J  and  substituting  in  (22)  we  shall  receive: 

a j 

3x  1/J  +  CCT 

4=0-' 


Fig.  1 1 . 


The  dependence  of  the  current  density  J  upon  time  t  at  B=(0;0,001;0,01;0,l;l;2)  is  shown  in 

As  we  see  availability  of  a  magnetic  field  strongly  influences  a  maximum  current  in  the  clot. 
Let’s  consider  now  the  solution  of  (22),  (22*)  at  T— >0.  Then: 


^o=1+l*Lt+K-1+fh 


0  T->0 


y  =  e(mp-l)Vat 

30  R 

—  =  pT 

3t 
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at 


J  =  e(mP-1^(l  -  a/x) 
e(/«p-i)t 

1  +  ax 

96  _  pe(mP~1)T 

3x-i  +  axe(w|3_1)T 
Then: 

If  1)  mb= 1 

6-1  =— ln(l  +  ax) 
a 


2)  mbi*l 


96 

9x 


>  0  always,  let's  consider  second  derivative 


926  pe^P-1^  ,  D  ,  N 

- T  =  7 - - - rr-(/«B-l-ax) 

^  (l  +  axe^P-^f 


The  sign  of  expression  in  brackets  determines  behavior  of  dependence. 


At  wP  —  1  >  a  dependence  represents  saturation. 

At  wP  - 1  <  a  there  is  an  inflection  point  at  Tper  = 


1  wP-1 
— — In— — 
wP-1  a 


It  is  seen  that  the  inflection  point  for  q  is  simultaneously  point  of  a  maximum  for  J. 
That  is: 
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J, 


max 


m(3  -  1/a 
l  +  ln(m[3-l/a) 


,  at  wP  - 1  >  a 


and 


*/max  =  1  (at  tr=0)  for  /Mp  - 1  <  a . 


Thus,  I\it  could  be  concluded  that  in  this  Section  a  simplified  model  of  clot  deformation  in 
turbulent  flow  has  been  developed.  Energy  balance  inside  the  clot  is  controlled  by  two  competing 
processes:  (1)  stretching  of  current  lines  due  to  exponential  growth  of  clot  surface,  and  (2) 
exponential  growth  of  conductivity.  Three  types  of  clot  temperature  growth  —  exponential,  linear, 
saturating  are  found.  Criteria  for  these  evolution  mode  has  been  obtained.  Limiting  growth  wall 
effect  was  shown.  The  main  turbulence  parameters  produced  by  magnetic  field  has  been  appritiated. 
Strong  influence  of  magnetic  field  to  maximum  clot  current  was  shown.  Results  of  proposed  model 
are  in  quantitative  agreement  with  experimental  observations. 

Current  Zero  at  Gaseous  Circuit  Breakers. 

Influence  of  Turbulence  [13]. 


It  can  be  easily  seen  that  the  behavior  of  the  electric  arc,  particularly  its  extinguishing  near 
the  zero  electric  current  point  at  the  alternative  voltage,  is  determined  of  course  by  the  kinetic 
processes  in  the  arc  column.  On  the  other  hand,  the  background  of  these  processes  and  according 
initial  conditions  are  to  a  considerable  extent  dominated  by  the  distributions  of  macroscopic 
parameters.  These  distributions  must  be  calculated  within  the  frames  of  hydrodynamic  description 
so  that  the  latter  is  a  very  important  part  of  the  problem  under  consideration.  Because  the  gas 
motion  in  the  considering  situation  is  turbulent  the  velocity  of  the  background  motion  is  a  random 
function  of  space  and  time  variables.  Unfortunately  the  intervals  of  time  typical  for  electrodynamic 
and  kinetic  processes  are  too  short  for  using  any  averaging  procedure,  so  the  standard  models  of 
fully  developed  turbulence  [19]  are  not  applicable  in  our  case  (neither  in  Euler  nor  in  Lagrange 
language  [15]).  The  alternative  way  in  which  some  «incomplete»  turbulence  can  be  taken  into 
account  was  proposed  in  [9]:  this  way  is  based  on  the  fact  that  the  material  lines  and  surfaces  at  the 
presence  of  the  random  velocity  field  grow  exponentially  (see  also  [15]).  Corresponding  model  of 
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the  exchange  of  energy  leads  to  the  equation  which  solution  shows  some  critical  behavior.  In  this 
connection  two  final  states  appear  to  be  possible  —  extinguishing  of  some  plasma  formation  (the  so 
called  clot)  or  its  developing  [10].  The  latter  shows  that  such  a  model  of  «incomplete»  turbulent 
phenomena  could  be  suitable  for  our  case,  of  course  after  generalization  of  the  model  with  regard  to 
the  nonhomogeneous  situation. 

In  addition  to  the  affirmations  mentioned  above  there  is  one  more  important  reason  for 
which  one  is  to  avoid  the  usual  procedures  of  the  statistical  or  time  averaging  (here  and  below  we 
denote  averaging  by  the  brackets  <  >).  Namely,  it  should  be  taken  into  account  that  electric 
conductivity  is  a  very  «strong»  (see  later)  function  of  the  temperature.  Consequently  one  can  write 
(at  the  nonuniform  temperature  distribution)  neither 

<«H<xr)> 

nor  (a)  =  (attj))}  , 

so  that  some  different  approach  for  the  finding  of  the  effective  values  of  the  conductivity  (or 
resistance)  must  be  found.  Because  the  conductivity  and,  consequently,  the  electric  arc  is  sharply 
non-uniform  the  main  problem  is  to  describe  the  behavior  of  such  strongly  non-uniform  plasma 
influenced  by  the  external  (non-stationary)  voltage.  This  problem  is  considered  as  follows. 

According  to  [9]  all  the  influence  of  the  turbulent  motion  on  the  arc  region  comes  to  taking 
into  account  of  the  exponential  growth  of  «material»  lines  and  surfaces  at  given  regions  of  the  arc. 
Such  a  viewpoint  on  the  random  velocity  field  influence  on  the  material  lines  and  surfaces  comes  to 
work  by  G.K.Bathelor  [14]  (this  work  included  an  error  which  was  remedied  later  in  [20], 
According  to  this  viewpoint  some  average  length  of  the  material  line,  for  instance,  depends  on  time 
in  the  next  form 

L(t)  =  L(t0)e^-^ 

where  E,  is  to  be  defined  as  a  value  of  order  of  the  «intrinsic»  frequency  of  turbulent  motion.  In 
[9,10]  slightly  different  quantity  was  used  which  depends  on  the  whole  range  of  turbulent 
frequencies  so  that  it  differs  from  E,  by  some  factor.  In  any  case  there  is  no  way  now  to  obtain  some 
«exact»  value  of  t,  in  both  approaches  and  as  a  final  result  it  is  to  be  estimated  from  the 


26 


experimental  data.  If  turbulent  motion  is  non-uniform  in  space  and  time  according  generalized 
relation  takes  the  obvious  form  (we  at  once  use  it  for  infinitesimal  element) 

dx  =>  <ixexp{jQ(x,t)^/}, 

which  will  be  used  below. 

Taking  account  of  the  fact  that  in  our  simulations  we  shall  watch  individual  fluid  particles  it 

is  easy  to  see  that  the  intrinsic  time  variable  t  is  to  be  given  in  the  form 

,  dx 

dt  - - 

u(x) 

where  w( x)  is  velocity  distribution  along  the  are  column.  This  time  variable  is  synchronized  with 
the  «extemal»  time  variable  (time  of  the  change  of  some  external  parameters,  for  instance,  voltage 
U (t) )  by  some  special  simulation  procedure. 

Returning  now  to  the  main  purpose  of  this  work  we  are  to  write  down  the  relations 
describing  the  mentioned  above  effects  of  turbulence  on  are  behavior. 

For  the  electric  resistance  of  the  arc  element  dx,  for  instance,  it  reads 

dR  <frexp{Q(x),t} 
a(x)5(x) 

where  a  is  the  electric  conductivity,  S  is  the  surface  of  the  electric  current  «tube»,  Q  is  the  time 
scale  factor  reflecting  the  influence  of  turbulent  motion;  its  dependence  on  distance  x  corresponds  to 
possible  non-uniformity  of  the  arc  region.  Namely  this  quantity,  £2(%),  is  to  be  connected  with 
turbulent  characteristics.  At  this  condition  the  Joule  dissipation  is 


q(x ,  t)  =  o(x)E2  0)exp{-  2Q(x)  •  , 

where  E(x)  is  (local)  electric  field.  Using  (*)  one  can  write  for  the  latter  (j  is  electric  current 
density) 


a 


U  (t)  exp{Q(x)  •  /} 


L 

o(x)Sq(x)  \dx 
0 


exp{2Q(x)  •  ?} 
o(x)S0(x) 
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where  the  symbol  «0»  corresponds  to  some  initial  conditions  concerning  current  distribution,  U ( t ) 
is  the  arc  voltage.  The  value  a  depends  explicitly  only  on  the  temperature,  T,  so  that 


a  =  Go  exP 


"  I 

f1  O" 

2kT 

U  T)_ 

where  I  is  the  ionization  potential,  k  is  the  Bolzmann  constant.  With  the  energy  equation  in  the  form 
dT 

pcp-^ 


and  the  equation  for  the  electric  circuit 

E(0  =  L—  +  IR 
dt 

where  L  is  the  inductance,  one  can  solve  the  problem  of  E,  T  and  j  definition,  in  the  column  of  arc. 
Introducing  some  characteristic  parameters  p*  (pressure),  /*,  R*,  assuming  L  to  be  the  length  of  the 
arc  and  coming  to  dimensionless  variables, 


q=q/q*  =  q  G0 


I*R 


*\2 


,  p  =  p/p*  ,  R*  =  -±-  ,  I*  =  Uq/R*  , 
O0S0 


-&  =  T/Tq  ,  G  =  g/gq  ,  m  =  l/2kT$  ,E-Ej 


I*  R* 


x  =  t/t*  (t*  =  L/R*)  ,  1  =  1/1*  ,  r  -R/R*  ,  Q  =  Q  t*  , 
x  =  x/L  ,  U  =  U/I*R*  ,  S=Sq/i}  ,  e=E/I*R*  , 


one  can  write  down  the  next  equations 


G  =  exp 

coin'd  _  y-1  t*  q 
dx  y  tq*p 


r  n 

m 

i — 

l  v\ 
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q  =  cE~  exp{-  2D  •  x}  , 
~  _  U  expjp  •  t} 

0  cSo 

e(x)  =  ^-  +  7r  , 
dx 


where  t^-p  *jq  * ,  y  is  specific  heat  ratio.  We  shall  start  with  the  solution  of  these  equation  in 
comparatively  simple  case  when  L=0.  At  this  condition  the  equation  (**)  takes  the  form 

d Inf}  _  7“1  q_ 

dx  y  P 

U(x)  =  7-r  , 


where  x  =  t/tg*  and  U(x)  is  some  given  function,  for  example,  U(t)  =  U$  sin  (Ot  where  coco  is  the 
frequency  of  the  alternative  voltage. 

Because,  on  one  hand  the  assumption  on  the  exponential  growth  of  the  material  lines  in 
turbulent  flow  renders  decisive  influence  on  the  arc  behavior  and,  the  other,  both  mentioned  above 
theories  of  this  effect  are  not  too  convincing  one  is  to  demonstrate  it  in  some  probable  way. 

Our  demonstration  is  based  on  the  representation  of  the  random  velocity  field  in  a  form  of 
sum  each  term  of  which  is  an  exact  solution  of  2-D  Euler  equation,  for  an  incompressible  fluid 

9v 

(W)V 

divV  =  0 
V  =  (u,  v) 

so  that  the  stream  function  \| is 

Vyt  =  Ak  cos {kx + ) cos (ky + yk)  , 

where  is  amplitude,  and  —  some  initial  values  of  the  phases.  It  is  easy  to  verity  that  \|/t 
satisfies  stationary  Euler  equation  in  the  form 
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By  9Ay  3\|/  3A\| f  _ 
dx  dy  dy  dx 

Consequently  we  construct  the  kinematics  of  turbulence  in  the  next  form 

dMf  n 

u  =  =  YukAk  cos(fcc + P* )sin(*v  +  Y* )  , 

^  k=\ 

v  =  — f  -  =  2  ■ Mt  sin(**  +  Pk )  cos(^  +  Yk )  • 
dx  *= i 

Placing  some  initial  «spot»  which  boundary  is  a  circle  which  diameter  D  >  2jt/&min  one 
can  observe  this  spot  developing  which  time  and  its  boundary  growth.  The  results  are  given  in 
Fig.  11.  From  the  latter  one  can  see  that  according  length  growth  is  really  exponential.  It  is  also  seen 
that  the  form  achieved  shows  its  randomness  except  some  regular  «cross»  (Fig.  11)  The  nature  of 
this  «cross»  can  be  clarified  as  follows.  Let  us  consider  the  set  of  the  point  at  which  u  =  v  (for 
instance,  we  consider  case  u,v  >0).  At  such  condition  the  relation 

sin(fcc + P^ )  =  cos  (kx + P^ ) 

and  similar  for  y  must  be  satisfied.  Therefore 

_  Jt 

kx  +  Qb  =  2nn  +  —  . 

™  4 
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Figll.  The  shapes  of  the  original  circular  contour  modified  by  a  ‘syntetic’  turbulence  field  with 
different  parameters. 


0.00  0.02  0.04  0.06 

Time,  ru. 


Fig.  12.  Confirmation  of  the  logarithmic  law  of  an  arbitrary  contour  length  growth. 


It  follows  from  the  last  relation  that  at  n  =  1  and  large  k  there  are  a  lot  of  points  with  u~v  in 
the  vicinity  of  the  origin  and  namely  these  points  give  rise  to  such  a  «coherent»  behavior  of  the 
boundary.  Next  we  intend  to  eliminate  such  behavior  by  introducing  some  temporal  factors  into 
velocity  definition  and  to  connect  the  rate  of  the  length  growth  to  «turbulent»  characteristics 
extracted  from  the  values  of  according  parameters. 
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It  should  be  noted  ones  more  that  qualitative  agreement  between  experimental  observations 
(see  Figs.2-4),  others  calculations  (see  Fig.5)  and  data  obtained  earlier  in  a  gas-plasma  flow 
(Figs.l  1)  confirms  the  similarity  of  the  processes  in  questions. 

The  most  important  practical  conclusion  is  that  so-called  synthetic  turbulence  can  be  used  in 
qualitative  analysis  of  miscellaneous  turbulence  effect  on  physical-chemical  processes  defined  by 
mixing  stage. 
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IV.  Reacting  Volume  Concept 


The  basic  idea  of  reacting  volume  concept  is  already  formulated  in  the  previous  Chapter. 

The  equation  of  energy  (6) 

P  Cp%  =  ^tV2T 

includes  the  heat  conductive  member  in  the  right  hand  side.  The  estimation  of  this  expression  used 
above  is  based  on  the  assumption  that  a  characteristic  value  of  the  temperature  gradient  VT  is 
defined  by  a  diffusion  length  //>.  The  latter  is  estimated  as  Id  ~y[Dt  where  D  and  t  diffusion 
coefficient  and  time,  correspondingly.  Time  t  is  calculated  from  the  moment  of  the  initial  contact 
between  the  hot  and  cold  fractions.  Because  of  well  known  physical  similarity  between  of  the 
conductive  heat  transfer  and  species  diffusion  it  is[possible  to  utilize  the  similar  approach  to  treat 
the  molecular  mixing  of  fuel  and  oxidant,  for  example,  by  means  of  diffusion  mechanism.  The 
generalized  driving  force  in  this  case  is  the  species  concentration  gradient  instead  of  the  temperature 
gradient  of  the  previous  case.  The  liquid  volume  formed  by  contact  surface  between  the  two 
gaseous  components  (fuel  and  oxidant,  for  example)  and  by  diffusion  thickness  Id  normal  to  this 
surface  growths  in  time  due  to  two  factors  at  least.  The  first  one  is  the  increasing  of  diffusion  length 

Id  ~  4~Dt .  In  non-uniform  non-steady  state  flows,  especially  in  turbulent  flows,  the  most  important 
is  the  second  factor  -  the  (exponential)  growth  of  the  contact  surface.  The  practical  model  of  such  a 
process  has  been  developed  in  the  study  of  a  current  carrying  plasma  clot  evolution  in  turbulent 
flows  and  illustrated  in  the  proceeding  examples.  Calling  such  an  intermediate  layer  as  reacting 
volume  one  can  attribute  to  the  fact  that  chemical  reaction  (burning  as  an  example)  between  two 
reactants  can  occur  typically  if  the  local  concentrations  both  reactants  are  high  enough.  Thus, 
reaction  can  start  only  inside  of  the  diffusion  layer.  Furthermore,  the  maximal  number  of  the 
reaction  acts  is  proportional  to  the  total  number  of  reactants  molecule  in  the  diffusion  layer.  For 
more  or  less  similar  gases  of  fuel  and  oxidant  the  upper  estimation  of  the  acts  is  about  one  half  of 
the  total  number  of  one  of  the  reactant  molecules  in  such  a  reacting  volume. 

It  is  clear  that  in  the  desirable  case  of  effective  burning  in  combustor  chamber  the  reacting 
volume  should  include  whole  amount  of  reactants  injected  into  the  chamber.  From  the  other  hand 
the  time  of  reacting  volume  formation  is  always  finite  and,  for  this  reason,  for  effective  combustion 
the  residence  time  should  be  longer  than  the  whole  reacting  volume  formation. 
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It  is  obvious  from  the  above  simplified  consideration  that  affecting  both  residence  time  and 
rate  of  reacting  volume  formation  can  result  in  combustion  efficiency  increasing. 

The  control  of  residence  time  is  out  of  the  scope  of  his  study.  Instead  the  methods  of 
intensification  of  the  reacting  volume  formation  are  a  main  goal  of  the  study.  Even  more  specific, 
the  new  ‘advanced’  technique  of  the  mixing  intensification  will  be  checked  -  this  is  the  MHD 
driven  intensification. 

The  background  of  this  technique  is  also  presented  in  the  previous  Chapter  when  speaking 
about  he  MHD  generation  (modification)  of  turbulence.  The  physical  background  is  the  non¬ 
potential  nature  of  the  jxB  body  force  that  directly  addresses  to  the  turbulence  (vorticity) 
production. 

The  additional  profit  of  he  MHD  technique  application  can  be  the  chemical  effects  due  to 
plasma  chemistry  mechanisms  which  can  take  place  because  of  rather  high  voltage  discharge  should 
be  used . 

In  more  detail  the  electrodynamics  and  maghnetoihydrodynamics  components  will  be 
described  in  the  next  Report. 


34 


V.  Concluding  Remarks  and  Future  Work 


1.  The  methods  of  kinematics  mixing  described  in  the  available  technical  papers  and 
textbooks  seems  to  be  rather  promising  for  the  aims  of  this  study. 

2.  The  experience  collected  in  the  non-uniform  gas-plasma  flows  study  ,  especially  connected 
to  the  plasma  clots  evolution  in  turbulent  flows  will  be  useful  for  the  high  speed  combustion 
studies. 

3.  The  estimation  of  the  expected  MHD  mixing  effects  has  shown  that  it  will  useful  for  the 
high  speed  combustion  applications. 

4.  Preliminary  analysis  of  the  experimental  validation  potentiality  resulted  in  positive 
conclusion. 

5.  The  adaptation  of  existing  numerical  models  for  turbulent  effects  treatment  in  gas-plasma 
application  to  the  combustion  case  is  successfully  started. 
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Introduction 


This  Report  #2  is  devoted  to  the  adaptation  of  the  existing  CFD  codes  and  chemical  kinetics 
code  to  the  case  of  interest  I  to  study  of  advanced  mixing  in  free  high-speed  jets.  The  preliminary 
plan  was  to  limit  this  study  only  by  mixing  study.  However  during  the  preliminary  analysis  it  was 
concluded  that  many  important  for  future  applications  phenomena  are  expected  under  combined 
kinematics,  dynamics  and  heat  releases  influence.  MHD  interaction  which  is  assumed  as  a  main 
newly  proposed  process  for  such  an  application  involves  (in  any  practical  system  based  on  thermal 
or  weakly  ionized  gases)  active  particle  such  as  ions,  electrons,  radicals,  exited  molecules  an  so  on 
extremely  important  for  the  high  temperature  combustion  especially  in  strongly  non-uniform 
systems.  Thus,  the  proper  chemical  kinetics  scheme  has  to  be  included  into  the  combined  numerical 
simulation  code.  For  these  proposed  we  decided  to  use  our  experience  and  numerical  codes  used 
earlier  to  study  the  processes  in  miscellaneous  plasma  chemistry  devices  and  in  internal 
combustion  engines.  The  KINEL  chemical  kinetics  code  was  developed  in  IVTAN  more  then  10 
years  ago,  checked  in  many  practical  applications.  This  code  has  been  also  modified  to  be  used  in 
the  combined  CFD  code  to  describe  multidimensional  chemical  reacting  gases  and  plasmas. 

Unfortunately,  the  full  scale  simulation  of  multi  dimensional  chemical  reacting  flow  with 
extremely  non-uniform  MHD  interaction  is  hardly  possible  in  nearest  future  with  the  computers 
available.  For  these  reasons  the  partial  simulation  is  planned  as  a  tool  to  understand  the  physical 
background  of  the  method  proposed.  Thus,  the  next  stage  will  be  a  simplified  physical  and 
corresponding  numerical  model  development.  The  approach  to  such  a  development  was  already 
formulated  in  previous  report  #  1,  where  he  concept  of  reacting  volume  was  formulated.  The 
reacting  volume  concept  has  to  proved  and  if  necessary  modified  with  the  data  obtained  with  more 
sophisticated  models. 

In  Section  1  the  general  physical  formulation  will  be  briefly  discussed. 

The  section  2  is  devoted  to  the  systematic  description  of  the  CFD  model  to  be  used  in  multi¬ 
dimensional  simulations. 

The  Section  3  presents  an  example  of  2  dimensional  time  dependent  simulation  of  the  electrical 
discharge  channel  moving  through  gaseous  media  due  to  body  force. 

The  section  4  presents  a  result  of  the  hydrogen/air  free  stream  combustion  with  a  simplified  9- 
componenets  species  scheme. 

In  the  last  section  future  steps  is  briefly  characterized. 
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1.  Mixing  at  Contact  Surface  in  Turbulent  Co-Flow  Streams 


As  an  example  of  practical  case  the  free  stream  jet  combustion  is  considered.  It  assumed  that 
two  parallel  jets  -  gaseous  fuel  and  air  of  equal  absolute  velocity  are  created.  The  main  point  of 
interest  is  the  contact  surface  between  these  two  components.  In  the  case  of  different  velocities 
the  Helmholz-Kelvin  instability  develops  very  fast.  Such  an  instability  process  provides 
turbulization  and  corresponding  turbulent  mixing  zone  development.  From  he  point  of  view  of 
free  stream  combustion  the  formation  of  turbulent  mixing  zone  is  extremely  positive  feature.  As 
we  already  discussed  in  the  previous  report  #1  the  turbulent  mixing  provides  in  fact  increasing 
the  effective  square  of  the  contact  surface.  Increased  contact  surface  results  in  intensification  of 
the  total  diffusion  fluxes  which  are  only  responsible  mechanism  to  prepare  flammable  fuel- 
oxidant  mixture  at  molecular  level. 

The  very  qualitative  vision  of  the  turbulent  mixing  development  at  the  contact  surface  is 
presented  in  Fig.  1  as  it  can  be  observed  in  moving  reference  frame.  The  mutual  penetration  one 
gas  component  into  the  other  one  is  characterized  in  time  by  multiplication  of  characteristics 
linear  scale.  It  was  shown  in  previous  study  that  such  a  scale  multiplication  results  in 
exponential  growth  of  he  contact  surface  in  particular. 


Fig.  1.  Schematic  of  ‘turbulent  mixing  ’  at  a  contact  surface  between  two  co-flow  streams 


However  despite  of  the  exponential  growth,  which  is  defined  by  a  characterization 
frequency  of  the  particular  turbulence  for  practical  application,  such  a  turbulence  mixing  can  be 
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not  enough  fast.  This  is  because  of  the  residence  time  is  very  often  very  short  for  high-speed 
combustion  chamber. 

To  intensify  the  contact  surface  growing  the  additional  mechanism  of  mixing  has  been 
proposed.  This  mechanism  involves  high  intensity  body  forces  acting  on  the  artificially  created 
electrical  current  channels  (arcs)  if  the  external  magnetic  field  of  a  proper  orientation  has  been 
applied.  As  it  will  be  demonstrated  in  the  Section  3  such  a  current  channel  can  be  shifted  across 
a  flow  from  one  component  (fuel  or  air)  to  the  neighboring  jet  (air  or  fuel,  correspondingly).  In 
small  time  scale  this  process  is  observed  as  “mushroom”  developing  across  the  contact  surface. 
From  the  point  of  view  of  the  multi  scale  turbulent  approach  MHD  powered  mixing  results  in 
modification  (in  the  proper  direction)  of  the  integral  characteristics  of  turbulence  such  as 
characteristic  turbulence  velocity  (frequency,  mixing  length  and  so  on). 

As  a  very  important  additional  advantage  of  the  MHD  forced  mixing  of  fuel-oxygen  system  is 
automatically  provided  ignition  and  combustion  stimulating  factors:  active  particles  (radicals) 
and  electrical  heating. 

In  order  to  check  the  real  effect  of  such  an  ‘advanced  mixing’  the  second  key  process  should  be 
explicitly  described  -  diffusion.  In  the  case  of  a  direct  simulation  of  a  single  mushroom 
development  it  can  be  provided  with  the  combined  model  described  in  the  next  section. 
However  for  multi  discharge  non-regular  case  we  plan  to  use  an  simplified  ‘engineering’  model 
based  on  several  semi-empirical  effective  parameters. 
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2.  Governing  equations 


Two-dimensional  Navier-Stokes  and  species  transport  equations  for  a  chemically  reacting 
gas  of  N  species  in  electric/magnetic  field  are  given  by 


3U  t  3F  |  3G  |  3FV  |  3GV  _£ 
3 1  dx  dy  dx  dy 


(1) 


where 


U  =  (p,pv;c)pv>,p£:,py,,L  ,pYNY 
F  =  lpvx,pv2  +  P,pvxvy,vx(pE  +  P),pvxYul  ,pvxYNJ 
G  =  {pvy,pvxvy,pv2y  +  P,Vy(pE+P),pvyYhL  ,pVyYNJ 
~Y^,Xxx’Xxy’Qx’J x  1>^*  * ^ xN / 

Gy  —  (fl>'Zxy’^yy’Qy’JyU^~  yN Y 

.  3vj.  2  „ 

Xxx  =  -2p.-^-  + — pVv 


Xxy  =  -V\ 


dx  3' 

<K+t>L 

dy  dx 


~  dvy  2 


avi+3v, 
dx  dy 


Vv  =  yix^^JL 


N  ,  \ 

6x  =  ~  vx^xx  vyXxy  ) 

dx  i= 1 

371  N  ,  \ 

Qy  ~  fyJyi  ~  (Vx^xv  ) 

9V  M 


^  =-PA«  ^  ,  i=U* 

OX 


(2) 
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In  these  equations,  the  physical  variables  are  the  density  p,  the  velocity  components  vx  and 

vy,  the  pressure  P,  the  temperature  T,  the  internal  energy  e,  the  total  energy  E  -e  +—{yl +  vy).  The 

mass-fractions  of  species  and  enthalpies  are  respectively  7,  and  /?,.  u  and  X  are  molecular  viscosity 
and  thermal  conductivity,  respectively,  of  the  gas  mixture.  The  diffusion  fluxes,  Jxl  and  Jyi  are 
defined  by  Fick’s  law  (2)  and  effective  diffusion  -  coefficients  of  species  are  obtained  as 

D  = — — ^ — 

I  *ilDv 
j*i 


where  ^ ^  WjWi  moje.fraction  of  i-th  species,  Wi  is  the  molecular  weight 

and  1=1  is  the  mixture  molecular  weight.  Binary-diffusivities,  Dy,  molecular  viscosity  and 
thermal  conductivity  are  evaluated  using  Chapman-Enskog  theory  and  Wilke’s  mixing  rule. 

Species  enthalpies,  hi,  are  expressed  as  follows 

hi(T)  =  hiJ+  J cpi(T)dT 

Tref 

where  as  the  formation  enthalpy  of  species  i,  cpi(T)  is  the  specific  heat.  Specific  heat  are 
approximated  by  5-th  order  polynomials  of  temperature. 

The  mixture  internal  enthalpy  and  energy  are  defined  as 


i=i 


e  =  h-  P/p  . 


The  pressure  is  determined  from  equation  of  stage  for  ideal-gas  mixture 


where  R  =  R()  jw ,  R°  is  the  universal  gas  constant. 

Source-term  accounts  contributions  from  chemical  reactions  for  species  transport,  forces 
from  electro-magnetic  interaction  and  the  work  of  electro-magnetic  field. 


S  =  (0,Fx,F  ,0£,w„— ,ww) 


Electro-magnetic  sources 

F  =  {Fx,Fy)=[jxB] 

9 

QE=(i-F)-Qrad  (3a) 


J  is  the  electric  current  density  and 
B  is  the  magnetic  induction  vector. 

E  is  the  electric  field  strength.  Ohm’s  law  is  used  to  link  E  and  J: 

J-i -[JxB]  =  c(E  +  [vxB])  (3b) 

B 

p  is  so  called  Hall  parameter  and  is  generally  approximated  by  function  of  pressure,  temperature 
and  magnetic  field.  Electrical  conductivity,  o,  can  be  approximated  as  function  of  temperature  of 

kind:  c(T) -  A-Tn  ■  e~v/r . 

For  the  slow  electro-magnetic  field,  (3b)  can  be  enclosed  by  definition  of  electric  scalar 
potential,  E=-grad  tp. 

Qrad  in  (3a)  denotes  the  radiative  losses  approximated  for  optically  thin  media  as 
Qrad  -  e  ■  aSB  '  ^  > 

(4)  where  £  is  the  effective  absorption  coefficient  and  Osis  is  Stefan-Boltzman  constant. 
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Chemical  sources 


The  species  rate  of  production,  w,  are  found  as 


if  any  of  Nr  reactions  can  be  written  in  the  form 

n  — N 

5>;,rx,  nr~  2X,*, 

i-i  *  i-i 

In  (3b)  v'  and  v"  are  the  stoichiometric  coefficients  of  forward  and  backward  reactions, 
respectively. 

c;  =  p  •  Yj  / W]  are  the  mole  concentrations, 

X]  =  W  ■  ci /p  are  the  mole-fractions. 

kfr  and  kbr  are  forward  and  backward  reactions  rate  constants,  respectively. 

Arrhenius  form  is  assumed  for  k/r  and  kbr- 

kjj.=a-Tb  ■  exp (-Ea/T) 

Solution  method 

System  of  equations  (1)  is  discretized  in  time  and  space  using  finite-volume  formulation: 

J  JU  da+  J(F  +  Fv)dSa  =  ]S  di o  (4) 

1  co  da 
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for  any  computational  cell  of  volume  to  enclosed  by  the  surface  S®  with  outward  unit  normal 
n(nx,ny). 


F  =  (F  +  Fv  )«*  +  (G  +  G  v 

If  piece-wise  constant  or  piece-wise  linear  distribution  of  variables  is  assumed  within  each 
control-volume  (cell),  the  transient  term  may  be  represented  as  co  •  ,  where  At  is  time-step. 

(4)  can  be  written  as  follows 


AU 

At 


+  R  =  0 


(4a) 


where  R  is  the  residual- vector. 

This  system  is  integrated  in  time  using  explicit  low-memory  Runge-Kutta  method: 


U0=U" 

U,-  =U0- a,  —  R(U,_,)  ,  I=m  (4b) 

CO 

u"*'  =  uK( 

where  U"  is  the  solution  on  time-level  f  and  U"+1  is  the  solution  on  time-level 


t"+1  =t"+At 

The  time-step,  At,  is  selected  as  minimum  of  convective  and  diffusive  time-scales  for  all  the 
mesh  control-volumes.  While  source-terms  for  momentum  and  energy  equations  are  computed 
explicitly,  calculation  of  chemical  species  production  rates  is  carried  out  implicitly  employing 
special  chemical-kinetics  solver. 


Integration  in  time  is  performed  by  the  program  KINEL  developed  in  IVTAN  Low- 
Temperature  Plasma  Division  for  solution  of  systems  of  ordinary  differential  equations  of  high 
stiffness  inherent  to  problems  with  chemical  and  thermal  non-equilibria.  On  every  time-step,  At, 
system  (4c)  is  integrated  based  on  “chemical”  time-steps  and  implicit  treatment  of  production-rates. 
The  chemical  time-steps  are  evaluated  within  KINEL  on  the  base  of  eigen-values  of  chemical 
source-term  Jacobian  matrix.  Low-cost  LU-decomposition  of  Jacobian  is  applied  on  each  chemical 
time-step  to  solve  the  system. 

Implicit  treatment  of  chemical  source-term  allows  one  to  avoid  strong  restrictions  to  the  gas- 
dynamic  time-steps  imposed  by  finite-rate  chemistry. 

The  integration  of  inviscid  and  viscous  flux-vector  over  the  control-volume  surface  is 
performed  as  follows.  One-integration-point  on  the  each  surface-element  (cell-interface  straddling 
two  neighboring  cells)  is  suitable  to  obtain  both  1-st  order  solution  (piece-wise  constant  variable 
distribution  within  the  cell)  and  2-nd  order  solution  (piece-wise  linear  variable  distribution)  in 
space.  The  inviscid  flux  at  the  integration-point  is  found  from  exact  solution  of  Riemann  problem 
(Godunov’  method),  i.e. 

¥ip  =  F(U;>),  and  Vip  =  Rie(VL,VR)  (4d) 

Rie(\Ji,  Ufl)  implys  solution  of  Riemann  problem  between  two  states  Ui  and  U* ,  left  and  right  state- 
vectors  on  both  sides  of  the  inter-face.  For  two  adjacent  cells,  P  and  E, 


U/,  -  Up  +\|/p(VUp,Arp^) 

U*=U£-\|f£(VU£,Ar£„)  (4e) 

VUp.p  are  the  cell-centered  gradients  and  ArPM  is  the  vector  directed  from  cell-center  P(E)  to  mid¬ 
point  of  the  interface  M.  the  gradient  in  each  cell,  VUp,  is  calculated  by  the  least-square  method  as 
minimum  of  functional 

-  ^L^jp  > 

J&(dp 

AJP  =  (Vj-lJP)-{VVp-ArPj) 
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Here  J  stands  for  any  cell  from  the  stencil  associated  with  cell  P.  usually,  five-point  stencil 
of  the  nearest  neighbors  (including  P )  is  considered.  The  limiter-functions  are  evaluated  to 
satisfy  the  monotonical  distributions  near  the  strong  flow  discontinuities.  In  the  smooth-flow 
regions,  the  second-order  accurate  spatial  distributions  are  recovered.  In  one-dimensional  problem, 
this  method  becomes  equivalent  to  the  well-known  Van  Leer  MUSCL  method. 

Viscous  fluxes  at  cell  interface  may  be  constructed  in  two  ways.  In  first  may,  the  face- 
gradient  is  taken  to  be  the  algebraic  average  ot  two  cell-gradients  described  just  above  with  adding 
up  the  projection  of  normal-derivative  to  prevent  the  non-physical  oscillations.  In  second  way,  the 
gradient  at  interface  is  calculated  by  application  of  Gauss’  theorem  to  auxiliary  cell” 


«>a(VUA/)=  J(U-nfl)dSa. 

9o)a 

Auxiliary  cell  represents  the  quadrilateral  vortices  of  which  are  two  cell-centers  (P  and  E) 
and  start-point  and  end-point  of  the  interface  segment  itself.  Both  approaches  were  found  to  be 
equally  appropriate. 

System  of  gouveming  equations  is  to  be  enclosed  by  the  initial  and  boundary  conditions. 
The  following  types  of  boundary  conditions  can  be  handled. 

1)  Inlet 

The  incoming  flow  compositions  is  specified  and  flow  variables,  p,  yJP  are  specified 
dependent  Mach  number  value.  For  supersonic  inflow  all  variables  can  be  set.  For  subsonic  inflow, 
total  temperature,  total  enthalpy  and  flow  angle  are  specified  and  iso-entropic  relations  are  assumed 
to  hold. 


2)  Outlet 

On  supersonic  outlet  all  variables  are  extrapolated  from  the  interior  with  1st  or  2nd  spatial 
accuracy.  Also,  static  pressure  can  be  specified  at  the  outlet  boundary. 

3)  Solid  Wall 
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Non-slip  conditions  are  applied  for  viscous  computations.  Adiabatic  or  constant-temperature 
conditions  are  used  for  energy  equation.  Non-catalytic  or  prescribed-compositions  conditions  are 
used  for  species  transport  equations. 

4)  Symmetry 

Normal-derivatives  are  set  to  zero  for  all  the  variables  except  normal  velocity  which  is  set  to 

zero. 


Depending  on  the  problem,  the  elliptic  problem  for  electric  potential,  which  is  the 
consequence  of  the  system  (3),  can  be  solved.  For  this  purpose,  high-performance  elliptic  solver  is 
employed  based  on  Additive  Correction  Multigrid  concept  along  with  the  Modified  Incomplete 
Choletski  methoid. 
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3.  Current  Channel  Propagation  in  External  Magnetic  Field 


The  evolution  of  the  initially  stated  uniform  arc  in  the  magnetic  field  was  computed  as 
demonstration  problem.  The  arc  was  assumed  to  be  uniform  in  z-direction  and  magnetic  field 
B=(0^y,0),  By-  1  Tesla,  was  specified  in  the  domain  6x1. 5mm  in  x-  and  y-directions,  respectively. 
The  uniform  arc  of  1mm  diameter  was  positioned  at  the  symmetry  axis  (y=0)  and  1.5mm  right  of 
the  left  boundary.  Initial  temperature  of  the  arc  has  been  set  to  9000K,  which  approximately 
corresponds  to  equilibrium  between  Joule  dissipation  and  radiative  cooling  in  the  uniform  arc. 

Initial  pressure  of  latm  has  been  specified  everywhere.  The  temperature  of  the  ambient  gas 
was  set  300K.  The  total  current  2A  was  assumed  to  flow  through  the  arc. 

From  these  considerations  the  two-dimensional  problem  on  evolution  of  the  arc  can  be 
solved.  Two  kinds  of  computations  were  carried  out:  inviscid  (and  thermally  non-conductive)  and 
viscous(thermally  conductive).  The  flow  fields  at  different  time  moments  are  shown  in  Fig.l 
(inviscid)  and  Fig.2  (viscous).  The  temperature  color-scale  ranges  from  300K  to  22000K  in  all 
figures.  Temperature  contours  are  uniformly  spaced  with  the  difference  440K. 

While  the  average  velocity  of  the  arc  is  nearly  same  for  both  calculations,  the  shape  of  it  is 
quite  different.  Inviscid  arc  constricts  in  its  leading  part,  whereas  trailing  part  is  strongly  cooled. 
Viscous  arc  is  more  uniform  during  time  period  considered.  However,  the  temperature  drops  below 


4000K. 


These  calculations  are  considered  as  the  preliminary  stage  of  the  whole  problem. 
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i 


R«tativ«  distance 

2.  Temperature  and  velocity  vectors  at  time-moments  10,20,30)40  and  50  mks  (from  up  to 
m).  Inviscid  case.  rmax=22000K,  A7’=440K. 


Relative  distance 


R*l atlv*  distance  Relative  distance  Relative  distance  Relative  distance 


4.  Reacting  Volume  Development  in  lDt  Simulation 


The  general  set  of  equations  was  reduced  to  lDt  case  to  simulate  a  simplified  case  of  reacting 
volume  development  in  an  axial  symmetric  configuration. 

As  a  typical  of scramjet  composition  the  example  of  hydrogen/air  system  is  considered. 

The  geometry  assumed  is  two  co-axial  configurations  consisting  of  internal  hydrogen  stream  and 
external  air  stream.  The  radius  of  internal  stream  is  1  cm,  and  outer  radius  of  external  air  stream  is  2 
cm. 

It  is  assumed  that  any  gas  parameter  is  changing  only  in  radial  direction  and  in  time. 

Initial  conditions  are:  uniform  temperature  T-  1200K,  pressure p  —  .  1  MPa  and  zero  velocity. 

The  boundary  conditions  at  the  outer  radius  r  =  2  cm  are  non-penetrated  wall  conditions. 

Simplified  widely  used  chemical  kinetics  scheme  involving  42  reaction  with  9  species  is 
implemented: 

1.  H2+02=0H+0H 

2.  0H+0H=H2+02 

3.  0H+H2=H20+H 

4.  H20+H=0H+H2 

5.  H+02=0H+0 

6.  OH+OH+02 

7.  0+H2=0H+H 

8.  0H+H=0+H2 

9.  0+H20=0H+0H 

10.  0H+0H=0+H20 

11.  H+H+M=H2+M 

12.  H2+M=H+H+M 

13.  H+02+M=H02+M 

14.  H02+M=H+02+M 

15.  0H+H+M=H20+M 

16.  H20+M=H+0H+M 

17.  H02+H2=H202+H 
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18.  H202+H=H02+H2 

19.  H02+H02=H202+02 

20.  H202+02=H02+H02 

21.  H+H02=0H+0H 

22.  0H+0H=H+H02 

23.  H+H02=H20+0 

24.  H20+0=H+H02 

25.  H+H02=H2+02 

26.  H2+02=H+H02 

27.  0+H02=0H+02 

28.  0H+02=0+H02 

29.  0H+H02=H20+02 

30.  H20+02=0H+H02 

31.  0H+0H+M=H202+M 

32.  H202+M=0H+0H+M 

33.  H02+H2=H20+0H 

34.  0H+H20=H02+H2 

35.  H02+H20=H202+0H 

36.  H202+0H=H02+H20 

37.  H+H202=H20+0H 

38.  H20+0H=H+H202 

39.  0H+M=0+H+M 

40.  O+H+MOH+M 

41.  0+0+M=02+M 

42.  02+M=0+0+M 

Under  conditions  considered  the  diffusive  combustion  develops.  The  typical  results  are  presented  in 
following  six  plots. 

In  Fig.4  the  radial  temperature  distributions  developing  in  time  are  presented.  It  is  clear  from  the 
plot  that  under  considered  conditions  the  combustion  is  diffusion  control  combustion  self  ignited  by 
the  initial  rather  high  temperature.  It  is  well  known  that  such  an  air-hydrogen  mixture  is  self-ignited 
when  the  temperature  higher  then  550-600C.  The  initial  temperature  assumed  in  this  example 
provides  high  combustion  rate.  The  maximal  temperature  about  3000  K  is  reached  in  .7  ms.  It  is 
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also  seen  that  at  the  late  stages  of  the  process  development  the  combustion  rate  decreases  due  to 
perhaps  the  increasing  of  diffusion  length. 


Fig.4  Development  of  Temperature  Distribution  across  Reacting  Volume 


The  heat  release  of  combustion  results  in  significant  disturbances  of  pressure.  At  he  initial  stage  of 
combustion  compression  waves  are  formed  at  the  reacting  volume.  The  Fig.5  illustrates  the  pressure 
distribution  development  in  time.  At  the  later  stages  the  pressure  changes  more  and  more  gradually 
corresponding  to  the  process  of  the  heat  release  in  a  closed  vessel. 
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Fig.5  Development  of  Pressure  Distribution  across  Reacting  Volume 


Fig.6  Development  of  Radial  Velocity  Distribution  across  Reacting  Volume 


The  pressure  disturbances  result  in  radial  velocity.  Because  of  the  boundary  conditions 

related  to  the  closed  volume  the  velocity  distribution  is  very  non-uniform  and  non-regular  in  time 
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(see  Fig.6).  It  should  be  noted  that  the  amplitude  of  the  radial  velocity  is  very  small  and  seems  to  be 
rather  insignificant  in  future  study  of  the  combustion  kinetics  in  free  stream  high-speed  jets. 


The  formation  of  the  reacting  volume  is  illustrated  in  Fig.7  where  the  time  variation  of  the 
radial  distribution  of  two  characteristics  species  are  plotted.  The  first  is  the  water  weight  fraction 
distribution,  and  the  second  -  the  nitrogen  weight  fraction  distribution.  The  first  value  is  a  main 
product  of  the  combustion  process  considered  here.  The  second  being  fully  non-changeable  due  to 
the  chemical  kinetics  scheme  used  represents  rather  accurate  the  diffusion  volume.  It  is  seen  that 
these  two  regions  are  very  closed  to  each  other  confirming  the  fact  that  in  this  particular  case  the 
combustion  is  the  diffusion-controlled  combustion.  A  slight  deformation  of  the  nitrogen  distribution 
at  the  outer  border  of  the  diffusion  volume  is  caused  by  he  the  weak  convection  and  characterizes 
the  secondary  effects  which  are  very  weak  here. 
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Fig.7.  Development  of  Radial  Distribution  N2  and  H20  across  Reacting  Volume 

The  full  species  information  is  presented  in  Fig.8.  The  radial  distributions  presented 
there  correspond  the  time  moment  of  500  mcs  when  the  combustion  riches  the  full  developed  stage. 
Fom  the  analysis  of  the  species  distributions  one  can  conclude  that  the  process  is  non-monotonic  in 
small  scale  consideration.  For  example  the  oxygen  distribution  reveals  an  intermediate  minimum 
due  to  combined  effect  of  diffusion  and  heating  by  heat  release.  At  the  same  time  we  recognize  that 
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such  a  non-notorious  behavior  is  observed  at  very  low  level  of  oxygen  concentration.  Under 
condition  different  from  the  considered  here  this  effects  can  be  much  more  significant. 


Fig.8  Main  Species  Distribution  across  Reacting  Volume  at  the  Time  Moment  of  500  mcs 

The  interesting  fact  is  presented  in  the  next  Fig.9  where  the  water  concentration  and  the 
temperature  are  presented  for  comparison.  The  non-trivial  results  here  is  that  the  peaks  of  the  water 
concentration  (combustion  product)  and  the  peak  of  the  temperature  resulting  from  heat  release 
from  the  combustiom  (water  formation)  are  shifted  of  each  other.  The  temperature  ‘propagates’  in 
front  of  water  formation.  Again  such  an  effect  can  be  significant  in  future  study  of  the  more  realistic 
case  of  combined  effects  of  gas  discharge  drowned  through  the  non-premixed  fuel/air  system. 
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Fig.9.  Comparison  of  Temperature  and  H2O  Weight  Fraction  Distributions  across  Reacting 
Volume  at  Two  Time  Moments  49  and  500  mcs. 

Concluding  this  example  presentation  it  could  be  recognized  that  the  chemical  kinetics 

simulation  prepared  for  the  future  implementation  into  combined  advanced  mixing  simulation  code 

satisfies  to  desirable  applications. 
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S.Concluding  Remarks  and  Future  Work 


Two  main  component  of  numerical  model  of  the  advanced  mixing  powered  by  MHD  body 
forces  are  developed. 

The  models  are  demonstrated  with  particular  cases  of  MHD  forced  mixing  and  with 
combustion  in  diffusion  mixing  area. 

The  next  step  is  a  verification  of  combined  model  involving  both  key  components. 

Basing  on  the  results  of  such  of  simulation  a  simplified  integral  model  of  the  MHD  forced 
mixing  and  combustion  in  reacting  volume  is  being  developed. 

It  was  also  found  even  on  this  preliminary  stage  that  the  MHD  mixing  has  provided  not  only 
intensification  of  the  mixing  (passive  combustion  enhancement)  but  also  the  dramatically 
increased  ignition  effectiveness  especially  for  practical  important  case  of  cool  fuel  jets. 
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Part  III 

1 .  MHD  Driven  Mixing  at  Contact  Surface  in  Co-Flow  Streams 

2.  Reacting  Volume  Development  in  2D  Simulation 

3.  Ignition  and  Combustion  under  MHD  Stimulation 

4.  Experimental  Validation  Scheme 

5.  General  Conclusion 
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1.  MHD  Driven  Mixing  at  Contact  Surface  in  Co-Flow  Streams 


As  a  typical  configuration  of  the  gas-  fuel  system  in  the  combustion  chamber  of  an  air-breathing 
engine  we  will  consider  two  co-flow  streams.  It  is  well  known  from  classical  fluid  dynamics  that  the 
contact  surface  between  two  co-flow  streams  is  absolutely  unstable.  The  Kelvin-Helmholz  instability 
destroys  the  contact  surface  resulting  in  multi  scale  chaotic  flow,  which  can  be  considered  as 
turbulent  flow.  It  seems  that  the  characteristic  size  of  perturbed  region  growths  downstream 
practically  linear  with  the  distance  from  the  initial  location  (nozzle  exit).  Moreover  one  can  assume 
that  the  rate  of  the  perturbed  region  cross-section  would  be  proportional  to  the  difference  of  the 
tangential  velocities  of  the  streams. 


Fig.l.  Schematic  of  MHD  driven  mixing  at  the  contact  surface  in  co-flow  streams. 

In  a  case  of  a  single  free  jet  we  consider  the  interface  surface  created  at  two  co-flow  streams  with  the 
tangential  velocities  difference  equal  to  the  free  jet  velocity  at  the  origin  cross-section.  For  such  a 
configuration  the  empirical  law  was  formulated.  This  law  says  that  characteristic  cross-section  size 
of  the  perturbed  region  growths  downstream  linearly  with  a  rate  factor  as  high  as  0.25.  We  can 
adopt  this  law  to  estimate  the  cross-section  characteristic  size  increasing  in  time  rather  than  in 

space.  Doing  so  we  can  find  that  the  rate  of  the  cross-section  size  growing  is  proportional  to  the  co- 
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flow  stream  relative  velocity.  This  rule  results  in  an  estimation  of  the  perturbed  region  characteristic 
cross-section  size  at  the  contact  surface  of  co-flow  streams  as  a  value  proportional  to  the  distance 
from  the  origin  location  with  the  factor  .25  X  whereof  is  a  reduced  relative  velocity 


v.  -  v. 


;r  = 


(v1+v2)* 


Thus,  the  maximal  rate  of  the  perturbed  region  growth  corresponds  to  the  free  jet  case,  and  minimal 
rate  corresponds  to  the  ‘neutral’  flow  when  the  velocities  of  both  streams  are  equal  to  each  other. 

It  should  be  noted  that  the  perturbed  region  mentioned  above  is  the  ‘turbulent’  mixing  region.  Only 
in  this  region  the  reacting  volume  can  develop. 

In  this  study  we  are  interested  to  assess  the  role  of  an  external  mechanism  stimulating  turbulence 
mixing.  For  this  purpose  it  is  convenient  to  consider  at  first  a  neutral  case  with  no  velocities 
difference.  Thus,  no  turbulence  mixing  will  occur.  Instead  only  MHD  driven  mixing  will  take  place. 
One  of  the  simplest  ways  to  provide  the  MHD  driven  mixing  is: 

a.  To  provide  electrical  discharge  channel  in  some  vicinity  of  the  contact  surface.  This 
discharge  channel  with  a  current  /  would  be  directed  along  the  contact  surface  and 
perpendicular  to  the  flow  velocity. 

b.  To  provide  external  magnetic  filed  of  the  magnetic  induction  B  co-oriented  along  the 
flow. 

c.  Thus,  the  body  force  jxB  normal  to  the  contact  surface  will  drive  the  current  channel 
across  the  flow.  Choosing  the  proper  current  direction  one  can  provide  the 
penetration  of  the  current  channel  from  the  one  stream  to  another.  Such  a  penetration 
results  in  extra  mixing  due  to  obvious  growth  of  the  area  of  the  contact  surface. 

Moreover  the  current  channel  created  and  supported  by  an  external  electric  field  is  an  region  with 
rather  high  energy  release  that  results  in  the  temperature  increasing  and,  if  the  electric  filed  high 
enough,  in  production  of  active  particles  -  radicals  which  is  primary  importance  to  initiate  the  chain 
combustion  reaction. 

Further  in  this  report  we  will  mine  the  first  two  results:  mixing  through  the  contact  surface  growth, 
and  ignition  of  combustion  due  to  heat  release  in  the  current  channel. 
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Fig.2.  The  mass  fractions  of  hydrogen  and  water  (left)  and  temperature(right)  for  Case  1  (no 
chemistry). 
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Fig.2.  The  mass  fractions  of  hydrogen  and  water(left)  and  temperature(right)  for  Case  1  (no 
chemistry,  cont.). 


In  the  Fig.2  of  this  section  the  reacting  volume  development  in  time  is  illustrated. 

In  configuration  of  the  Fig.  1  with  no  relative  velocity  the  problem  of  the  MHD  driven  mixing  is 
considered.  It  is  assume  that  total  current  of  a  single  discharge  channel  as  high  as  10  Amps;  the 
external  magnetic  field  is  3  T.  the  gas  properties  are  assumed  to  be  correspondent  to  the  air- 
hydrogen  flow.  In  this  particular  case  no  chemical  reaction  is  considered.  At  the  same  time  all 
diffusion,  heat  conductance  and  viscosity  effects  are  included  as  it  has  been  defined  in  Part  II 
numerical  model.  The  initial  static  pressure  is  1  atm.  The  initial  gas  temperature  is  300K,  the  initial 
temperature  inside  the  arc  -  discharge  channel  is  assumed  to  be  equal  to  6000K  that  more  or  less 
characteristic  temperature  for  the  considered  conditions.  The  arc  channel  diameter  is  2  mm,  the 
initial  location  of  the  discharge  channel  in  2  mm  from  the  contact  surface.  In  the  figs  hydrogen,  and 
the  upper  occupy  the  lower  part  of  the  calculated  domain  by  air. 

The  time  period  of  the  first  70  mcs  has  been  considered.  It  can  be  easily  seen  from  the  pictures  that 
the  diffusion  starts  immediately  creating  a  reacting  volume.  The  rough  estimation  of  such  a  process 
made  in  Part  I  of  this  report  is  qualitatively  confirmed.  In  time  the  convection  influence  becomes 
more  and  more  significant.  The  square  of  the  interface  surface  increases  especially  when  the 
discharge  current  has  penetrated  the  interface.  In  case  of  discharge  initiation  in  the  hydrogen  (very 
light  gas)  the  discharge  channel  deforms  significantly  while  approaching  the  heavy  gas  border 
spreading  along  the  interface. 

Probably  due  to  rather  low  velocity  of  the  discharge  channel  the  diffusion  region  at  the  tip  of  the 
cross-section  of  the  current  duct  differs  not  significantly  from  those  at  the  peripheral  region  of  the 
domain.  It  means  that  the  reacting  volume  can  be  estimated  as  a  product  of  the  interface  surface 
square  Sc  and  diffusion  thickness  Sn~vbt.  More  accurate  conclusion  at  this  stage  of  study  is 
problematic  because  of  many  mechanisms  affecting  the  formation  of  the  reacting  volume  are 
involved.  For  example,  the  temperature  dependence  of  the  diffusion,  real  gas  properties,  and  others 
factor  are  obviously  important. 

As  it  was  expected  the  nature  of  the  MHD  enhancement  of  the  kinematics  mixing  is  the  vorticity  of 
the  Lorentz  body  force.  The  flowfield  presented  in  the  previous  Part  II  is  clearly  intensively  vortex. 
The  vorticity  (simplified  turbulence)  results  in  fast  increasing  of  any  fluid  surface  located  in  such  a 
flow. 
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The  heat  release  effect  in  an  applied  electric  field  is  illustrated  by  the  second  column  in  Fig.2  where 
the  temperature  filed  evolution  is  presented.  It  can  be  observed  that  the  area  of  the  hot  gas  increases 
as  well. 

In  the  next  section  we  will  examine  the  reaction  effects  on  the  mixing  and  combustion  under 
conditions  considered. 
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The  Fig.3  represents  the  process  of  MHD  driven  mixing  and  ignition.  The  result  has  been  obtained 
with  full  numerical  model  described  in  the  previous  Part  II.  The  model  includes  2Dt  simulation  of  the 
reacting  gases  with  accounting  of  diffusion,  heat  conductance  and  viscosity  effects.  The  chemical 
kinetics  scheme  used  has  been  also  described  in  the  previous  Part.  This  scheme  describes  the 
combustion  of  hydrogen  in  air  on  the  base  of  a  simplified  scheme  consisting  of  42  reactions  with  nine 
species  involved.  Nitrogen  is  considered  as  a  passive  component  considered  as  a  third  body  molecule. 
At  this  stage  no  ions  chemistry  has  been  included.  Thus,  the  results  obtained  in  this  simulation  should 
be  considered  as  qualitative. 

In  the  left  column  of  Fig.3  the  weight  fraction  of  the  hydrogen  and  water  -  main  final  product  of  the 
combustion  are  presented.  The  hydrogen  concentration  level  lines  (dark  monocolor  lines)  indicated  the 
diffusion  region  where  the  hydrogen  fraction  differs  from  1  (the  initial  value  in  a  lower  part  of  the 
domain).  Due  o  rather  low  temperature  in  gas  surrounding  the  arc  the  rate  of  chemical  reaction  is 
negligible  low,  no  combustion  is  observed  in  primary  diffusion  region.  In  the  contrast  to  this  the 
diffusion  region  formed  at  the  tip  of  the  arc  cross-section  moving  to  the  contract  surface  is  heated-up 
by  the  energy  release  in  the  arc.  So  the  temperature  of  the  mixture  here  is  high  enough  to  initiate  the 
hydrogen  oxidation  reactions.  The  water  fraction  is  represented  in  Fi.3  by  colored  level  lines.  It  is 
clearly  seen  that  since  the  arc  approached  the  interface  surface  he  water  production  increases  rapidly. 

It  is  interesting  fact  that  maximum  water  concentration  does  no  coincide  with  the  maximal  temperature 
location  but  always  behind  of  that  from  the  hydrogen  side. 

Future  study  will  be  devoted  to  parametric  investigation  of  such  a  process  of  MHD  driven  ignition. 
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4.  Experimental  Validation  Scheme 

Schematic  Diagram  of  Experimental  Setup 
to  demonstrate  MHD  Advanced  Mixing 
in  Co-Flow  Free  Streams 


w 


+  I  1 

Fig.4  Schematic  diagram  of  an  experimental  setup 
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In  Fig.4  one  of  the  possible  arrangements  of  small  laboratory  setup  to  demonstrate  the  effect  of  MHD 
driven  mixing  and  combustion  in  co-flow  streams  is  presented.  The  facility  includes  two  gas  sources 
(‘fuel’  and  ‘oxidant’),  magnetic  induction  system,  and  electrical  discharge  system.  The  magnetic  filed 
induced  along  the  main  velocity  of  jets  in  combination  of  the  arc  discharge  organized  near  the  contact 
surface  of  the  two  jets  create  the  body  force  driving  the  arc  channel  across  the  free  jets.  According  the 
preliminary  estimates  and  numerical  simulation  such  a  configuration  ca  improve  the  combustion 
performance  significantly. 

In  more  detailed  the  proposed  experimental  facility  will  be  described  in  the  Final  report  of  the  related 
project  supported  by  ISTC. 
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General  Conclusion 

The  new  concept  of  MHD  driven  Advanced  mixing  of  gaseous  fuel-oxidant  system  in  a  high-speed 
combuster  are  formulated. 

The  background  of  such  a  concept  consists  of  the  earlier  developed  model  of  kinematics  mixing  in 
vortex  turbulent  flows. 

The  vortex  nature  of  the  Lorentz  body  force  resulting  from  perpendicular  vectors  of  electrical  current 
and  externally  applied  magnetic  filed  provides  the  flow  vorticity  intensification. 

The  enhancement  of  the  mixing  under  such  conditions  is  characterized  by  an  exponential  growth  of 
the  contact  surface  square  in  a  turbulent  flow. 

The  second  mechanism  -  molecular  diffusion  is  responsible  for  the  creating  the  flammable  mixture  in 
an  physically  infinitely  small  volume. 

These  two  mixing  mechanisms  are  combined  to  formulate  the  concept  of  reacting  volume  -  the 
volume  where  only  the  combustion  is  possible. 

Besides  of  the  qualitative  physical  model  the  proper  numerical  model  has  been  developed  to  study  the 
process  on  the  microscopic  scale. 

This  model  includes  the  treatment  of  chemical  reacting  flows  on  the  base  of  frill  Navier-Stockes 
approach,  chemical  kinetics  model,  and  MHD  effects  model. 

Preliminary  results  confirm  the  basic  physical  ideas.  The  combustion  in  co-flow  streams  can  be 
effectively  initiated  and  control  by  the  manipulation  of  he  electrical  parameters:  strength  of  external 
magnetic  field  and  electrical  discharge  current. 

The  preliminary  sketch  of  the  future  laboratory  scale  experiment  is  presented  as  well. 

The  parametric  study  of  the  proposed  concept  as  applied  to  the  conditions  typical  for  he  practical 
situation  in  a  scramjet/  ramjet  will  be  completed  in  a  related  project  supported  by  ISTC  in  2000/2001. 
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